﻿#include "Self_Define_Functions.h"
#include <time.h>
#include <cstdio>
#include <algorithm>
#include <iostream>
#include <fstream>
#include <stdlib.h>
#include <iomanip>
#include <string>
#include<cstring>
#include <boost/random.hpp>
#include <boost/random/uniform_int.hpp>
#include <boost/random/cauchy_distribution.hpp>
#include <boost/random/normal_distribution.hpp>
#include <boost/random/uniform_real.hpp>
#include <stdio.h>
#include <math.h>
#include <malloc.h>
#define PI 3.1415926535897932384626433832795029
double *OShift, *M, *y, *z, *x_bound;
int ini_flag = 0, n_flag, func_flag, *SS;
using namespace std;

int main(int argc, char *argv[])
{
	//////////////////////////global processing/////////////////////////////////////

	int funToRun[] = { 1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30 };  //function set

	int funNum = 30; //total function num
	double optima[30] = { 100,200,300,400,500,600,700,800,900,1000,1100,1200,1300,1400,1500,1600,1700,1800,1900,2000,2100,2200,2300,2400,2500,2600,2700,2800,2900,3000 };

	/////////////////////////////////////////////////////////////////////////////////////

	int i, j;
	double derta1, derta2;
	double p;
	double pp1 = 0.15;
	double dd1 = 0.1;
	int NP = 150;
	NewType *sorted_data = new NewType[NP];





	double p_J = 0.2;//The parameter of JADE top best individual
	int top_gbest_num = ceil(NP * p_J);
	int A_num = 0;//The number of individuals in Archive
	double Erfa = 0.5;
	double p_num; //Elite group
	double *c_strategy_pro = new double[8];//Probability of each mutation strategy being selected according to convergence
	double *d_strategy_pro = new double[8];//Probability of each mutation strategy being selected according to diversity
	double *select_mutation_pro_c = new double[8];//Save the probability for roulette wheel selection 
	double *select_mutation_pro_d = new double[8];//Save the probability for roulette wheel selection 
	double sum_con;
	double sum_dis;
	int *select_mutation = new int[NP];//record the  selected mutation strategy
	double global_fitness_old;
	double mid_CR;
	double mid_F;
	double CR;
	double F;
	double MIN = -100;//Low bound
	double MAX = 100;//Up bound
	double *CR_record = new double[NP];
	double *F_record = new double[NP];
	double *global_best = new double[dim];
	double global_fitness;
	int *sorted_index = new int[NP];
	int *individual_rank = new int[NP];
	double **population = new double*[NP];
	double **A = new double*[NP];//Archive
	for (i = 0; i < NP; i++)
	{
		population[i] = new double[dim];
		A[i] = new double[dim];
	}


	double *results = new double[NP];// the fitness results for the whole population
	double *trail_vector_results = new double[NP];

	int FV = 0;
	int function_index;
	int run_index;
	int parent1, parent2, parent3, parent4, parent5;

	double *mutation_vector = new double[dim];
	double **trail_vector = new double*[NP];
	for (i = 0; i < NP; i++)
		trail_vector[i] = new double[dim];

	int print_cont = 1;//the number of record best fitness
	int FV_ind = MAX_FV / 10;//record interval
	for (function_index = 0; function_index < funNum; function_index++)
	{
		cout << "Function " << funToRun[function_index] << " Begined!" << endl;

		boost::mt19937 generator(time(0)*rand());

		// to initialize the population
		boost::uniform_real<> uniform_real_generate_x(MIN, MAX);
		boost::variate_generator< boost::mt19937&, boost::uniform_real<> > random_real_num_x(generator, uniform_real_generate_x);

		boost::uniform_real<> uniform_real_generate_r(0, 1);
		boost::variate_generator< boost::mt19937&, boost::uniform_real<> > random_real_num_r(generator, uniform_real_generate_r);

		//select an individual from the population
		boost::uniform_int<> int_generator1(0, NP - 1);
		boost::variate_generator< boost::mt19937&, boost::uniform_int<> > int_number1(generator, int_generator1);
		boost::uniform_int<> int_generator2(0, top_gbest_num);
		boost::variate_generator< boost::mt19937&, boost::uniform_int<> > int_number2(generator, int_generator2);

		char fun1[10];
		char fun2[10];
		char fun3[10];
		snprintf(fun1, 10, "%d", funToRun[function_index]);
		snprintf(fun3, 10, "%d", dim);
		snprintf(fun2, 10, "%d", NP);
		string filename_fitness = "./N" + string(fun2) + "_dim" + string(fun3) + "/Fitness_result_for_function_" + string(fun1) + ".txt";
		ofstream out_fitness(filename_fitness.c_str());

		if (!out_fitness)
		{
			cerr << "Can not open the file " << filename_fitness << endl;
			exit(1);
		}

		for (run_index = 0; run_index < timesOfRun; run_index++)
		{
			cout << "Running the " << run_index + 1 << "th times" << endl;
			FV = 0;
			A_num = 0;
			print_cont = 1;
			//Initialization
			for (i = 0; i < NP; ++i)
			{
				for (j = 0; j < dim; j++)
				{
					population[i][j] = random_real_num_x();
				}

				cec14_test_func(population[i], &results[i], dim, 1, funToRun[function_index]);
				results[i] = results[i] - optima[function_index];

				FV++;
			}
			Sort_Population_Fitness(sorted_index, individual_rank, results, NP, sorted_data);
			global_fitness = results[sorted_index[0]];


			while (FV < MAX_FV)
			{
				//Calculate the parameter for mutation select
				derta1 = 0.1 + dd1 * pow((1.0*FV / MAX_FV), Erfa);
				derta2 = 1 - derta1;
				p = 0.05 + pp1 * pow((1.0*FV / MAX_FV), Erfa);
				p_num = p * NP;
				sum_con = 0;
				sum_dis = 0;
				for (i = 0; i < 8; i++)
				{
					//calculate the weight for each muation strategy 
					c_strategy_pro[i] = 1 / (sqrt(2 * PI)*derta1 * 8)*exp(pow(7 - i, 2) / (2 * derta1*derta1 * 8 * 8)*(-1));
					d_strategy_pro[i] = 1 / (sqrt(2 * PI)*derta2 * 8)*exp(pow(i, 2) / (2 * derta2*derta2 * 8 * 8)*(-1));
					sum_con += c_strategy_pro[i];
					sum_dis += d_strategy_pro[i];
				}
				//calculate the probability for each muation strategy 
				for (i = 0; i < 8; i++)
				{
					c_strategy_pro[i] = c_strategy_pro[i] / sum_con;
					d_strategy_pro[i] = d_strategy_pro[i] / sum_dis;
				}
				//calculate the probability for roulette wheel selection 
				select_mutation_pro_c[0] = c_strategy_pro[0];
				select_mutation_pro_d[0] = d_strategy_pro[0];
				for (i = 1; i < 8; i++)
				{
					select_mutation_pro_c[i] = select_mutation_pro_c[i - 1] + c_strategy_pro[i];
					select_mutation_pro_d[i] = select_mutation_pro_d[i - 1] + d_strategy_pro[i];
				}

				Sort_Population_Fitness(sorted_index, individual_rank, results, NP, sorted_data);
				global_fitness_old = global_fitness;
				global_fitness = results[sorted_index[0]];
				memcpy(global_best, population[sorted_index[0]], sizeof(double)*dim);
				for (i = 0; i < NP; i++)
				{
					CR = (individual_rank[i] + 1.0) / NP;
					F = (individual_rank[i] + 1.0) / NP;
					boost::mt19937 generator(time(0)*rand());
					boost::normal_distribution<> Cr_generator(CR, 0.1);
					boost::variate_generator< boost::mt19937&, boost::normal_distribution<> > CR_number(generator, Cr_generator);
					boost::normal_distribution<> F_generator(F, 0.1);
					boost::variate_generator< boost::mt19937&, boost::normal_distribution<> > F_number(generator, F_generator);
					mid_CR = CR_number();
					mid_F = F_number();
					while (mid_CR <= 0 || mid_CR >= 1)
					{
						mid_CR = CR_number();
					}
					CR_record[i] = mid_CR;
					while (mid_F <= 0)
					{
						mid_F = F_number();
					}
					if (mid_F >= 1)
						mid_F = 1;
					F_record[i] = mid_F;
				}

				for (i = 0; i < NP; ++i)
				{
					//Sort_Population_Fitness(sorted_index, individual_rank, results, NP,sorted_data);
					boost::mt19937 generator1(time(0)*rand());//A_num在进化开始过程中值是改变的，所以这是选取P并A中的随机个体的随机数
					boost::uniform_int<> int_generator4(0, NP + A_num - 1);
					boost::variate_generator< boost::mt19937&, boost::uniform_int<> > int_number4(generator1, int_generator4);
					if (individual_rank[i] < p_num)
					{
						select_mutation[i] = Roulette_Selection(select_mutation_pro_c, 8) + 1;
					}
					else
					{
						select_mutation[i] = Roulette_Selection(select_mutation_pro_d, 8) + 1;
					}

					switch (select_mutation[i])
					{

					case 1:
						parent1 = int_number1();
						while (parent1 == i)
						{
							parent1 = int_number1();
						}

						parent2 = int_number1();
						while (parent2 == parent1 || parent2 == i)
						{
							parent2 = int_number1();
						}

						parent3 = int_number1();
						while (parent3 == i || parent3 == parent1 || parent3 == parent2)
						{
							parent3 = int_number1();
						}
						parent4 = int_number1();
						while (parent4 == i || parent4 == parent1 || parent4 == parent2 || parent4 == parent3)
						{
							parent4 = int_number1();
						}
						parent5 = int_number1();
						while (parent5 == i || parent5 == parent1 || parent5 == parent2 || parent5 == parent3 || parent5 == parent4)
						{
							parent5 = int_number1();
						}
						Mutation_rand_2(population[i], mutation_vector, population[parent1], population[parent2], population[parent3], population[parent4], population[parent5], MIN, MAX, F_record[parent1], dim);
						break;

					case 2:
						parent1 = int_number1();
						while (parent1 == i)
						{
							parent1 = int_number1();
						}

						parent2 = int_number1();
						while (parent2 == parent1 || parent2 == i)
						{
							parent2 = int_number1();
						}
						parent3 = int_number1();
						while (parent3 == parent2 || parent3 == parent1 || parent3 == i)
						{
							parent3 = int_number1();
						}
						Mutation_rand_1(population[i], mutation_vector, population[parent1], population[parent2], population[parent3], MIN, MAX, F_record[parent1], dim);
						break;
					case 3:
						parent1 = int_number1();
						while (parent1 == i)
						{
							parent1 = int_number1();
						}

						parent2 = int_number1();
						while (parent2 == parent1 || parent2 == i)
						{
							parent2 = int_number1();
						}

						parent3 = int_number1();
						while (parent3 == i || parent3 == parent1 || parent3 == parent2)
						{
							parent3 = int_number1();
						}
						parent4 = int_number1();
						while (parent4 == i || parent4 == parent1 || parent4 == parent2 || parent4 == parent3)
						{
							parent4 = int_number1();
						}
						parent5 = int_number1();
						while (parent5 == i || parent5 == parent1 || parent5 == parent2 || parent5 == parent3 || parent5 == parent4)
						{
							parent5 = int_number1();
						}
						Mutation_current_to_rand_2(population[i], mutation_vector, population[i], population[parent1], population[parent2], population[parent3], population[parent4], population[parent5], MIN, MAX, F_record[i], dim);
						break;

					case 4:
						parent1 = int_number1();
						while (parent1 == i)
						{
							parent1 = int_number1();
						}

						parent2 = int_number1();
						while (parent2 == parent1 || parent2 == i)
						{
							parent2 = int_number1();
						}
						parent3 = int_number1();
						while (parent3 == i || parent3 == parent1 || parent3 == parent2)
						{
							parent3 = int_number1();
						}
						Mutation_current_to_rand_1(population[i], mutation_vector, population[i], population[parent1], population[parent2], population[parent3], MIN, MAX, F_record[i], dim);
						break;

					case 5:
						parent5 = int_number2();
						while (sorted_index[parent5] == i)
						{
							parent5 = int_number2();
						}
						parent1 = int_number1();
						while (parent1 == i || parent1 == sorted_index[parent5])
						{
							parent1 = int_number1();
						}
						parent2 = int_number1();
						while (parent2 == parent1 || parent2 == i || parent2 == sorted_index[parent5])
						{
							parent2 = int_number1();
						}
						parent3 = int_number1();
						while (parent3 == i || parent3 == parent1 || parent3 == parent2 || parent3 == sorted_index[parent5])
						{
							parent3 = int_number1();
						}
						parent4 = int_number4();
						while (parent4 == i || parent4 == parent1 || parent4 == parent2 || parent4 == parent3 || parent4 == sorted_index[parent5])
						{
							parent4 = int_number4();
						}
						if (parent4 < NP)
						{
							//size of Archive <= NP
							Mutation_cur_pbest_2(population[i], mutation_vector, population[i], population[sorted_index[parent5]], population[parent1], population[parent2], population[parent3], population[parent4], MIN, MAX, F_record[i], dim);
						}
						else
						{
							//size of Archive > NP
							Mutation_cur_pbest_2(population[i], mutation_vector, population[i], population[sorted_index[parent5]], population[parent1], population[parent2], population[parent3], A[parent4 - NP], MIN, MAX, F_record[i], dim);
						}
						break;

					case 6:
						parent3 = int_number2();
						while (sorted_index[parent3] == i)
						{
							parent3 = int_number2();
						}
						parent1 = int_number1();
						while (parent1 == i || parent1 == sorted_index[parent3])
						{
							parent1 = int_number1();
						}
						parent2 = int_number4();//parent3为Xr2(有存档)
						while (parent2 == i || parent2 == parent1 || parent2 == sorted_index[parent3])
						{
							parent2 = int_number4();
						}
						if (parent2 < NP)//分两种情况，A没有NP个个体时与A已经有NP个个体时
						{
							Mutation_cur_pbest_1(population[i], mutation_vector, population[i], population[sorted_index[parent3]], population[parent1], population[parent2], MIN, MAX, F_record[i], dim);
						}
						else
						{
							Mutation_cur_pbest_1(population[i], mutation_vector, population[i], population[sorted_index[parent3]], population[parent1], A[parent2 - NP], MIN, MAX, F_record[i], dim);
						}
						break;

					case 7:
						parent1 = int_number1();
						while (parent1 == i || parent1 == sorted_index[0])
						{
							parent1 = int_number1();
						}
						parent2 = int_number1();
						while (parent2 == parent1 || parent2 == i || parent2 == sorted_index[0])
						{
							parent2 = int_number1();
						}
						parent3 = int_number1();
						while (parent3 == i || parent3 == parent1 || parent3 == parent2 || parent3 == sorted_index[0])
						{
							parent3 = int_number1();
						}
						parent4 = int_number1();
						while (parent4 == i || parent4 == parent1 || parent4 == parent2 || parent4 == parent3 || parent4 == sorted_index[0])
						{
							parent4 = int_number1();
						}
						Mutation_current_to_best_2(population[i], mutation_vector, population[i], global_best, population[parent1], population[parent2], population[parent3], population[parent4], MIN, MAX, F_record[i], dim);
						break;

					case 8:
						parent1 = int_number1();
						while (parent1 == i || parent1 == sorted_index[0])
						{
							parent1 = int_number1();
						}
						parent2 = int_number1();
						while (parent2 == i || parent2 == parent1 || parent2 == sorted_index[0])
						{
							parent2 = int_number1();
						}
						Mutation_current_to_best_1(population[i], mutation_vector, population[i], global_best, population[parent1], population[parent2], MIN, MAX, F_record[i], dim);
						break;
					}
					Crossover(trail_vector[i], population[i], mutation_vector, CR_record[i], dim);
					cec14_test_func(trail_vector[i], &trail_vector_results[i], dim, 1, funToRun[function_index]);
					trail_vector_results[i] = trail_vector_results[i] - optima[function_index];
					FV++;
					if (trail_vector_results[i] < results[i])
					{
						if (A_num < NP)//|A|<NP
						{
							memcpy(A[A_num], population[i], sizeof(double)*dim);
							A_num++;
						}
						else//|A|==NP
						{
							memcpy(A[int_number1()], population[i], sizeof(double)*dim);
						}
						results[i] = trail_vector_results[i];
						memcpy(population[i], trail_vector[i], sizeof(double)*dim);
						if (results[i] < global_fitness)
						{

							global_fitness = results[i];
							memcpy(global_best, population[i], sizeof(double)*dim);
						}
					}
				}
				if (FV > FV_ind * print_cont)
				{

					out_fitness << global_fitness << endl;
					print_cont++;
				}
			}

		}
		out_fitness.close();

		cout << "Function " << funToRun[function_index] << " Finished!" << endl;
	}


	//release the resouces
	for (i = 0; i < NP; ++i)
	{
		delete[]population[i];
		delete[]trail_vector[i];
		delete[]A[i];
	}
	delete[]select_mutation_pro_c;
	delete[]select_mutation_pro_d;
	delete[]c_strategy_pro;
	delete[]d_strategy_pro;
	delete[]population;
	delete[]A;
	delete[]trail_vector;
	delete[]individual_rank;
	delete[]results;
	delete[]trail_vector_results;
	delete[]mutation_vector;
	delete[]CR_record;
	delete[]F_record;
	delete[]global_best;
	delete[]sorted_index;

	return 0;
}



